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Abstract 

We study the application of the density matrix renormalization group 
(DMRG) to systems with one-dimensional acoustic phonons. We show how 
the use of a local oscillator basis circumvents the difficulties with the long- 
range interactions generated in real space using the normal phonon basis. 
When applied to a harmonic atomic chain, we find excellent agreement with 
the exact solution even when using a modest number of oscillator and block 
states (a few times ten) . We discuss the use of this algorithm in more complex 
cases and point out its value when other techniques are deficient. 
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In a previous paper |l|], we have studied how the real-space density matrix renormaliza- 
tion group (DMRG) || could be applied to dispersionless phonons in the ground state 
of magnetic or fermion chains. We showed that the unavoidable truncation of the local 
oscillator space proved manageable for a Peierls or spin-Peierls deformation. Acceptable 
numerical accuracy could be obtained by properly selecting the dimension of the oscillator 
space as a function of the amplitude of the lattice deformation (the Peierls gap). 

The true Peierls, spin-Peierls, and superconductivity problems, or more generally the 
electron-phonon problem, however, involve acoustic phonons having dispersion. Are the 
hydrodynamic modes problematic? Let us look at other numerical approaches and see 
how they deal with phonons. Although the dispersion has not created any difficulty in 
exact diagonalization calculations on a one- dimensional (ID) lattice [|J, most efforts with 
electron-phonon systems have used dispersionless phonons || or have focused on a single 
acoustic mode related to an order parameter ||. The reason was to keep the dimension 
of the quantum space tractable. The chain or cluster size being rather small (of order 10 
sites), however, there is little hope of doing numerical justice to the hydrodynamic modes. 
The quantum Monte Carlo method (QMC) J7| has apparently done very well with acoustic 
phonons and the hydrodynamic modes. The QMC does have sign problems, however, with 
the electronic part in frustrated systems and with some magnetic impurity problems |§. 
There is no problem with phonons, of course, in a momentum space formulation of the 
DMRG ||. The treatment of local interactions (Hubbard interaction, exchange interaction, 
electron-phonon interaction...) is, however, more complex and computationally exacting in 
momentum space. 

The DMRG has the capability of reaching chain lengths in the hundreds, is free of the sign 
problem of the QMC, and is ideally suited to short-range interactions. The situation with 
the hydrodynamic modes is potentially troublesome, however. The long-range displacement- 
displacement autocorrelation function ((ue — u^) 2 ) in an infinite chain increases logarithmi- 
cally with the distance Rfgi between the static (classical) atomic equilibrium positions. The 
same is true of the displacement amplitude squared ((ue) 2 ) at the center of a finite chain 
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with closed boundary conditions (fixed ends) which grows logarithmically with the number 
N a of atoms in the chain. As a consequence, the amplitude of the atomic displacements can 
grow without bounds as the chain length increases in the DMRG procedure. In other words, 
there is no fixed point for the infinite-chain DMRG algorithm. It is forever fleeting as the 
chain grows. There is also another constraint on the form the Hamiltonian should have in 
the DMRG. The procedure takes place in real space and, ideally, the Hamiltonian should be 
reducible to a form involving only short-range coupling terms between the central site(s) and 
the two adjacent blocks. Otherwise, the numerical accuracy and memory expenditures - for 
the storage of the block variables coupling to the central site(s) - suffer greatly. It so turns 
out, as we shall see below, that using the normal phonon coordinates, defined in momentum 
space, leads to long-range couplings in real space and thus to potential difficulties for the 
DMRG. 

It is the purpose of this paper to study the constraints imposed by these considerations 
on the DMRG treatment of the ground state of a chain bearing acoustic phonons. We have 
chosen to focus specifically on the phonons leaving aside any electronic or spin counterpart. 
There are well documented applications of the DMRG to spin || [TIJ [[□]], electron [T2| 



or mixed [ITSf systems. As such, the atomic chain problem is more a testing ground for 



the DMRG than a physical problem since its excitations and thermodynamics have exact 
solutions. 

We shall study finite chains having N a atoms and closed boundary conditions. Keeping 
the end atoms at fixed positions eliminates the troublesome q = mode which, for periodic or 
open boundary conditions, corresponds to uniform displacements (and unbounded {{ui} 2 )). 
We shall assume the usual harmonic Hamiltonian 



and the boundary conditions u\ = un s =0. The solutions are of course well known in 
terms of the annihilation d q and creation d^ q operators for the normal modes of momentum 
q = rnrN~l d , 1 < n < N moc i, where N mo d = (N a — 2) is the number of normal modes. One 



has simply 




g 



where the eigenfrequency is u q 



2^y K/m |sin(g/2)|. Here and throughout, we have put the 



static equilibrium interatomic distance equal to one. 

The need for writing the Hamiltonian in real space stems from the real space algorithm 
used in the DMRG. In the case of our Eq. |2|, one could try to exploit the operators d,£ = 



The F(x) function is slowly decreasing, varying as x~ 2 , at long distances. Consequently, the 
block operators Ou 1 — F(l' — £)de' are relevant for all sites £' in the blocks. The storage 
requirements are large. The number of stored elements is of order N^M^ at each step s, for 
the infinite-chain algorithm, where Nt, s and are the number of atoms and the number of 
selected states in each block, respectively. It is huge, of order ^2 s N a M£, for the finite-chain 
algorithm, the sum covering all steps s of the procedure. The numerical accuracy would 
greatly suffer. 

There is a second problem, for the infinite-chain algorithm, related to the quantization 
of the phonon momentum q. Its values depend on the chain length N a and thus change 
from step to step. Consequently, so does the amplitude F(x) of the coupling terms. This 
is rather annoying. This problem can of course be circumvented by using the finite-chain 
algorithm as one can use the q values of the chosen chain length. 

The way around these problems is to use another basis set of quantum oscillators. Let 
us use instead the local oscillators (as in ||) that are solutions of 



N m od J2 g Gxp(iq£) d q on the perhaps natural reflex of using the diagonal phonon annihilation 
operators. Unfortunately, the Hamiltonian in real space then has the form 




(3) 



F ( X ) = N rnld Yl ^1 eX P(^) • 
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where u a = y/2K/m, 



2mw r 



(b\ + be). With this construction, one can write 



H = H Q + H c 



(5) 



where 



UiUe+i 
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£ 5>l + fc)(&J*i + fc+i) • 



This expression for H is exact. The coupling terms H c are now, however, short ranged and 
independent of the chain length. This last formulation is thus highly preferable for use with 
the DMRG even though the local oscillators are not diagonal. As a matter of fact, the 



connection between the b e and the d q is not trivial: 



h = Nj d J2 [(w /2u; ff )5 + (u q /2u )*\ sm(q£)d q 
q 



(6) 



This is a canonical transformation. We shall now study the implementation of Eq. [|. 

As mentioned in the Introduction, the central atom(s) oscillator space must be truncated 
at a dimension sufficient to properly represent the vibrational motion of the atoms. Regard- 
less of the algorithm used, one must begin by using the infinite-chain algorithm which is 
plagued by the absence of any fixed point. Indeed, the vibrational amplitude of the central 
atom was already mentioned to be logarithmically increasing at each step: 

<M 2 > = N-l d J2 (J^) (Mdi)) 2 cx ln(iV a ) . (7) 

This last quantity is a measure of the average occupation number (n^) of the local oscillator. 
It is then a prerequisite that the dimension M u of the local oscillator space be sufficient to 
cover the requirements of Eq. 0, that is M v 3> (ne). Our previous experience [J] has revealed 
that this dimension should be of order ten. We have thus used only one central atom (site) in 
order to keep the computation time and storage requirements within acceptable values. We 
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have monitored the occupation probability P s {n) of each oscillator state |n) = (n\)~^(b\) n \Q) : 
< n < M u , of the central atom for each step s of the infinite-chain algorithm for iV a up to 
49 atoms. We find that P s {n) ~ 0.35 exp(— n/a s ), a s ~ .33 + .181n(iV a — 3), fits the data 
quite well for most occupations, except for the first and last ones n = 0, M v — 1 . Notice 
the logarithmic dependence for a s which weighs the average and the standard deviation 
for the local oscillator occupation. Having this information, one can estimate a minimal 
relative error on the oscillator space statistics 5P S ~ P S (M U — 1) at each step resulting from 
using a finite dimensional oscillator space for the central atom. This error is propagative 
and thus, the total minimum numerical error of a sweep should be expected to be of order 
8P « P S (M U ). One can thus use this criterion to determine the minimal requirement on 
M v for the situation at hand. 

As the finite-chain method leads to the greatest precision let us first look at some 
results using this algorithm. Table | lists a few runs for a 25 atom chain. We have found 
that two sweeps are sufficient in all situations. Further sweeping changes little to the nu- 
merical values. Table | shows the various parameters used in each run: the oscillator space 
dimension M u , the number of block states kept M^, the number of target states M t , and the 
estimated minimal error 5P arising from the truncated oscillator space. We have calculated 
the numerical error on the ground state correlation energy 5E corr = (E® Q ^ RG — E° orr ) / E° orr , 
where E corr = E gs — \N mo dfvjJ and E gs is the ground state energy. The superscripts 
"DMRG" and "o" refer to the computational and the exact values, respectively. Note 
that E° s = | Y2q tWq- This correlation energy is more significant than E gs since the zero 
point energy of the local oscillators is not a numerically meaningful quantity. We have also 
calculated the error on the correlation function between the central atom and its first neigh- 
bor 5C = \((u t - u^i) 2 ) DMRG - ((ui - u e _i) 2 ) o \ I ((u e - u £ _i) 2 ) o , the one on the oscillation 
amplitude square of the central atom 5U = \{( u i) 2 ) dmrg ~ (M 2 )ol / (M 2 )o> and the er- 
ror on the average oscillator occupancy of the central atom 5n = \(n) o — (n) DMRG \ / (n) . 
(n) a = \b\btj , while ug and bg are defined in Eq. |] and || respectively. The run with 11 
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oscillator states and Mb = 10 shows SE corr to be roughly of the same size as SP. Increasing 
My does not change SE corr much. A saturation level has been reached for M v > 11 even 
though 5P keeps on decreasing. The central oscillator space is sufficiently large that en- 
larging it any further does little. Why is this? The reason is that important information 
contained in the central site is not being relayed to the next iteration step because too few 
block states are kept. Indeed, doubling Mb, as in the last run, considerably reduces the error 
on the ground state energy. One notices that the error on the nearest-neighbor correlation 
is of the same order as 5E corr . The errors SU and Sn, however, are of the same order but 
much larger than SE corr . This is systematically found in all our simulations. The reason 
is that these quantities sample the local oscillator excitations to a much larger degree than 
the other quantities. The accuracy on excitation energies is less than on the energy in the 
DMRG. Similar findings, shown in Table [H], were obtained for 49 and 99 atom chains. 



In Table |T|, we used the infinite-chain algorithm to look at the effect of increasing the 
number of target states. We see that the accuracy gets better when using a few target states 
and worsens when using too many. This is quite typical of the DMRG [TTJ. Using four 
target states transfers better information, through the projection procedure, to the block 
states and improves 5n considerably. Using too many target states dilutes the relevant 
information with irrelevant one. The runs with 30 block states show that the infinite-chain 
algorithm with four target states can perform as well as the finite-chain one (compare to the 
last entry in Table |) although the latter has just a single target state and a smaller number 
of block states. One then has to weigh the numerical expenditure of having more states as 
compared to doing two sweeps. This result is most important as it shows how to make a 
proper use the infinite-chain algorithm. 

We have shown that expressing the Hamiltonian in terms of local oscillators leads to a 
short-range coupling between blocks and central atom(s). This is a nice property to have 
when implementing the DMRG procedure. We have observed that the numerical error 
introduced by the truncation of the local oscillator space can be satisfactorily controlled by 
ab-initio selection of the size of the oscillator space for the problem at hand. This applies to 



any chain-length envisaged. We have found that the number of block states to be used 
is quite critical to the accuracy. As a rule of thumb, one should keep at least something like 
2 times the dimension of the local oscillator space M v . More important though, we have 
shown that the infinite chain algorithm can be used in spite of the anticipated limitations 
with regard to the lack of a fixed point and the logarithmic dependence of the atomic motion 
as a function of chain length. 

We have also found the usual DMRG truncation error 8E trunc is useless as an error 
indicator. We monitored it and observed it to be much smaller than the error on the 
correlation energy 5E corr or even 5P. This is rather puzzling since, as seen in the previous 
section, important errors stem from the quality of the projection of the target states onto 
the block states. The reason has to do with the limited information contained in the ground 
state of the super block with regard to the central site. The central site of the super block 
is an inversion symmetry center whereas it no longer is so at the end of the new block. 
The symmetry breaking information that is required for a proper description of the last site 
of the new block is thus contained in block states that have little statistical weight in the 
density matrix. 

The Hamiltonian we have studied is of an academic nature. By choosing the local 
oscillator basis, we knowingly sacrificed its phonon conservation law in order to obtain a 
form that is in harmony with the DMRG procedure. Our study does, however, provide a 
necessary stepping stone for problems that do not conserve the number of phonons. These 
are the more physical situations like the Peierls or spin-Peierls chains or, quite generally, 
the electron-phonon problem in one dimension. The computational expenditure required 
for phonons is rather modest and would presumably remain so in electron-phonon or spin- 
phonon problems with a gap as in [[L]]. 

We acknowledge the financial support of the Natural Sciences and Engineering Research 
Council of Canada and the Fonds pour la Formation de Chercheurs et d'Aide a la Recherche 
of the Quebec government. 
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TABLES 

TABLE I. Finite-chain algorithm for 25 atoms. The parameters are: the oscillator space 
dimension M u , the number of block states the number of target states Mt, the estimated 
minimal error on the oscillator statistics 8P, the numerical error on the ground state correlation 
energy one on the correlation function between the central atom and its first neighbor 

5C, the one on the oscillation amplitude squared of the central atom 8U, and the error on the 
average oscillator occupancy of the central atom Sn. 
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TABLE II. Finite-chain algorithm for chains of ^=49 and 99 atoms. The parameters are 
defined in Table I. 
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TABLE III. Infinite chain algorithm for 25 atoms. The parameters are defined in Table I. 

M v M b M t SE corr 5n 

15 20 1 3.7 x 10" 4 0.0035 

15 20 4 2.8 x 10" 4 0.0015 

15 20 7 1.7 x 10~ 3 0.0085 

15 30 4 2.6 x 10" 5 3.5 x 10" 4 

15 30 1 1.5 x 10" 4 .00145 
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